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ABSTRACT 

Lower  bounds  are  obtained  for  the  covariances  of  the 
orbital  elements  of  an  initially  circular  orbit  due  to  omissions 
of  harmonica  of  a  given  degree  in  the  expansion  of  the  earth's 
gravity  potential.  This  is  accomplished  by  applying  the 
methods  of  linear  estimation  theory  to  the  surface  gravimetry 
results  of  Kaula,  extended  to  satellite  altitudes.  Numerical 
results  are  presented  for  the  case  of  a  100-n  mi 
geocentric  orbit. 
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SECTION  I 


INTRODUCTION 


With  the  increased  accuracy  requirements  placed  upon  the 
investigator  in  the  field  of  astrodynamics,  it  has  become  necessary  to 
r ~ tcir.  ^r><»ff^cj«nt8  in  the  expansion  of  the  earth's  gravity 

potential  in  terms  of  spherical  harmonics.  However,  the  inaccuracies 
inherent  in  omitting  harmonics  of  a  given  degree  have  yet  to  be  estimated. 
This  analysis  provides  lower  bounds  for  these  effects  for  the  class  of 
initially  circular  geocentric  orbits. 


a 
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THE  COVARIANCE  EQUATIONS 


It  is  assumed  that,  for  small  deviations  from  circularity,  the  elements 
of  an  initially  circular  geocentric  orbit  are  governed  by  a  linear  system  of 
differential  equations  of  the  form 

d  =  F(t)  AX  (t)  +  G(t)u(t)  (i-a) 


where 

X  =  a  6  X  1  column  vector  representing  the  two  body  orbital 
elements 

u(t)  =  a  3  X  1  column  vector  of  zero  mean  representing  the 
forcing  function 

F(t)  and  G{tj  =6x6  and  6X3  matrices,  respectively 
A  =  a  deviation  from  circular  valuep 
t  =  the  time 

The  solution  to  this  system  of  equations  is 

t 

t)  G(t)  u(t)  dT  (i-b) 

o 


AX(t)  a  *(t, 


t0}AX(tc) 


♦J 


where 

I 

■V 

o  =  a  subscript  notation  indicating  an  epoch  value 

$>(t,  t)  =  the  6  6  state  transition  matrix,  which  transforms  perturba¬ 

tions  in  the  orbital  elements  at  time  r  into  perturbations  at 
time  t. 


Now,  since  this  analysis  is  limited  to  initially  circular  orbits,  we  have  from 
the  definition  of  A 


AX(tQ)  =  0 


Thus,  Eq.  (1-b)  can  be  rewritten 


AX(t)  = 


t)  G(t)  u(t)  dT 


Further,  noting  that  u(t)  is  of  zero  mean,  from  Eq.  (1-c) 


(1-c) 


COV[AX(t)]  =  r  f  ^(t,  t)  G(t)  E  |u(t)  UT(q)jGT(q)  0T(t,  rj)  dT  dtj  (2) 


where 

E  =  an  expected  value 
COV  =  the  covariance 

The  quantity  X  is  now  equated  to  the  following  set  of  orbital  elements 
(see  Fig.  1): 

p  =  radial  separation  between  the  satellite  and  geocenter 
v  -  speed  of  the  satellite 

P  =  angle  between  the  radius  vector  to  the  satellite  and  the  velocity 
vector 

i  -  the  inclination 

fl  -  longitude  of  the  ascending  node 

£  =  angle  in  the  orbit  plane  between  the  point  of  maximum  declination 
and  the  instantaneous  satellite  location  (positive  in  the  direction  of 
satellite  motion)  i.  e.  ,  the  argument  of  latitude  minus  t/2  radians. 
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equatorial 

PLANE 


Figure  1.  Orbital  Elements 


It  is  convenient  to  replace  time  by  £,  as  the  independent  variable  (and 
conversely  to  replace  4  by  t  as  an  orbital  element)  thus  modifying  Eq.  (Z)  to 

COV[AX(U]  =  j  J  <t>(L  t)  G(t)  E  [u(t>  uT(r1)j  GT(q)  ^T(t.  t,)  dt  dn 
^o  ’’o 

(3) 


Prior  to  direct  integration,  expressions  for  the  terms  appearing  in  Eq.  (3)  as 
functions  of  the  orbital  elements  and  the  components  of  the  forcing  function 
must  be  found.  The  analysis  to  accomplish  this  purpose  follows. 


SECTION  in 


EVALUATION  OF  TERMS  APPEARING 
IN  THE  COVARIANCE  INTEGRALS 


Neglecting  second  and  higher  order  terms  in  AX  yields  in  the  absence 
of  the  forcing  function  {Ref.  1) 


Extending  Eq.  (4)  to  include  the  effects  of  the  perturbing  accelerations  as 
given  by  the  Gaussian  equations  of  motion  with  £  replacing  time  as  independent 
variable  (see  e.  g.  ,  Ref.  2)  yields 


-  7- 


=  v  a  p  +  —  s 

dt  a  H  Vo 


dAi 

dC 


p  sin 
ro 


W 


(5) 


aas  -!°f!!Lw 

dt  2  , 


v  sin  i 
o 


d  At 


Pn  PD  MB  C  cos  i 

Ap  -  — ~ -  Av  + - j - W 

v  v  sin  i 

o  o 


where  R,  S,  and  W  are  the  components  of  the  perturbative  acceleration;  R 
lying  in  the  direction  of  the  increasing  radius  vector,  W  in  the  direction  of  the 
moment  of  momentum  vector,  and  S  in  the  instantaneous  orbit  plane  with  sense 
such  that  R,  S,  and  W  form  an  orthogonal  right-handed  set.  It  should  be  noted 
that  Eq.  (5)  is  equivalent  to  the  matrix  Eq.  (i-a)  with  £  replacing  t  as  the  inde¬ 
pendent  variable  and 


AX 


Ap 

Av 

AP 

Ai 

A» 

At 


(6-a) 


8- 


Since  all  of  the  coefficients  appearing  in  Eq.  (4)  [or  Eq.  (6-b)]  are  constants, 
it  is  a  straightforward  matter  to  solve  this  Bystem  of  equations  for  the  transi¬ 
tion  matrix  yielding 


f 

Z  -  cos(t  -  t) 

?  P0  r  i 

[l  -  co.(f,  -  t)] 

o 

-P0  -  t) 

0 

0 

<r 

r2,  [co«(l  -  t)  -  1] 

P0 

2  cos  -  t)  -  1 

vo  «tn(t  -  t) 

0 

0 

0 

-i-  .ln(l  -  t) 
po 

-  •£•  .in(t  -  t) 
a 

co«(t  -  t) 

0 

0 

0 

>-  (6-e) 

0 

0 

0 

1 

0 

0 

0 

0 

0 

0 

1 

0 

~C3(t  -  r)  -  Z  »ln(C  -  t)] 

[j(t  -  T)  -  4  »in(t  -  r)] 

2Po  r 

-  - —  D  -  coi(;  -  t)] 

0 

0 

The  R,  S,  and  W  components  of  the  correlation  matrix  of  the  forcing 
function,  due  to  effects  of  harmonics  of  given  degree,  are  derived  in  Appendix 
A.  This  is  accomplished  by  extending  the  surface  gravimetry  results  of 
Kaula  (Ref.  3)  to  satellite  altitudes.  Examination  of  these  results  shows  all 
three  components  to  be  linear  combinations  of  cosines  of  integer  multiples 
of  the  satellite  separation  angle.  Further,  assuming  these  components  are 
not  cross-correlated  yields  a  correlation  matrix  with  terms  of  the  form 
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0 


E  [u(t)  uT(t)0  s= 


ot  cos  X.j(t  -  tj) 
0 

G 


0 

P  cos  \2(t  -i\) 
0 


0 

"Y  COB  X^(t  “  T 


where  a,  p,  and  the  integers  Xj.  X.,.  X,  are  defined  in  Appendix  A 


SECTION  IV 


EVALUATION  OF  THE  COVARIANCE  INTEGRALS 


It  should  be  observed  that  a  difficulty  associated  with  the  limits  placed 
upon  Eq.  (3)  is  the  implicit  assumption  that  the  covariance  of  the  orbital 
elements  is  correlated  with  a  monotonically  increasing  separation  angle  t,. 

Of  course,  separation  central  angles  are  never  greater  than  it  radians.  For 
pu -poses  of  the  present  analysis,  this  difficulty  is  ameliorated  by  neglecting 
correlations  beyond  a  preselected  separation  angle  denoted  by  A  (where  A  <  tr). 
This  necessitates  reformulating  the  right  hand  side  of  Eq.  (3)  as  the  three 
double  integrals 


COV  Ux((|] 


-t,  +A  .  . 

j° 


I(t,  q)  dT  dq 


"A  rT+A 

+  f  I  I(t,  q)  dT  dq 


J  / 

JC+A  ■'t  -i 


(7) 


+ 


I(t,  q)  dT  dq 


where 

I(t,  q)  =  The  integrand  of  Eq.  (3)  =  ^(t,  v)  G(t)  E |u(t)  uT(q)J  GT(q)  £T(t,  q) 

The  (shaded)  area  in  the  t,  q  plane  over  which  Eq.  (7)  is  integrated  is  shown  in 
Fig.  H.  Of  course,  Eq.  (7)  holds  only  when  £  -  2  2A.  It  should  also  be 


noted  that  a*  formulated  Eq.  (7)  yields  only  a  lower  bound  for  the  covariance 
of  the  orbital  elements  since  it  completely  neglects  the  covariance  arising 
from  re-crossings  of  previous  ground  traces,  (1.  e.  ,  it  is  equivalent  to  the 
assumption  that  the  satellite  passes  over  new  topography  with  each  revolu¬ 
tion).  Substituting  Eqs.  (6-d),  (6-e),  and  (6-f)  into  Eq.  (7)  yields,  after 
direct  integration,  the  symmetrical  matrix 


COV 


Ap 

all 

a12 

a13 

0 

0 

a16 

Av 

a21 

a22 

a23 

0 

0 

a26 

Ap 

a3i 

a32 

a33 

0 

0 

a36 

Ai 

0 

0 

0 

a44 

a45 

a46 

Afl 

0 

0 

0 

a54 

a55 

a56 

At 

. 

_a6t 

a62 

a63 

a64 

a65 

a66^ 

(8) 


where  the  matrix  elements  are  given  in  Appendix  B. 


SECTION  V 


NUMERICAL  EXAMPLE 


In  order  to  perform  the  summation  in  Eqs.  (A- 14),  (A- 17),  and  (A- 18), 
numerical  values  are  needed  for  the  degree  variances  tr£  .  For  this  purpose, 
the  results  obtained  by  Kaula  in  Refs.  (3)  and  (4),  for  the  degree  variances  of 
the  terrestrial  free -air  anomalous  gravity  field  degree  variance,  are 
adopted.  These  values  are  illustrated  as  a  function  of  degree  n  in  Fig.  3. 
Unfortunately,  Kaula  has  only  determined  degree  variances  through  n  =  32. 

At  n  =  30  he  obtained  a  negative  degree  variance,  a  result  which  is  theoret¬ 
ically  impossible  as  indicated  by  Eq.  (A- 13).  Hence,  for  this  analysis,  degree 
variances  beyond  n  =  30  were  represented  by  =  35.  09  (0.  949)n  (a  rela¬ 
tionship  chosen  so  as  to  approximate  degree  variances  below,  and  be  well 
behaved  above,  n  =  30).  The  solid  line  in  Fig.  3  illustrates  degree  variances 
obtained  from  this  relationship.  Below  n  -  30,  Kaula's  degree  variances 
were  utilized.  Further,  correlations  beyond  separation  angles  of  it  radians 
were  neglected,  i.  e. ,  A  was  set  equal  to  tt. 

The  degree  variances  described  above  were  then  utilized  to  obtain 
numerical  results  for  the  diagonal  terms  of  the  matrix  comprising  the  right 
hand  side  of  Eq.  (8).  In  Figs.  4  through  9,  the  graphs  of  the  square  root 
of  these  diagonal  elements  (representing  the  standard  deviations  of  the 
orbital  elements)  are  presented  as  a  function  of  revolution  for  the  case  of  a 
satellite  in  a  circular,  100-n  mi  polar  orbit.  A  small  oscillatory  component 
superimposed  on  the  secuh.r  trend  has  been  smoothed  out  of  these  figures. 

Each  of  these  figures  consists  of  three  curves  corresponding  to  three 
assumptions  regarding  the  gravity  model.  The  curve  with  the  largest  ordi¬ 
nates  was  obtained  by  summing  from  n  =  2  to  n  =  50;  that  with  the  second 
largest  ordinates  by  summing  from  n  =  7  to  n  =  50;  and  that  with  the  smallest 
ordinates  by  summing  from  n  =  13  to  n  =  50.  The  exclusion  of  the  contribu¬ 
tion  of  the  lower  degrees  terms  in  a  spherical  harmonic  representation  of  the 
geopotential  is  equivalent  to  the  assumption  of  a  perfect  knowledge  of  the 


KAULA  S  VALUES  REPRESENTED  BY  + 
EXTRAPOLATED  VALUES  BY  SOLIO  LINE 
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Figure  3.  Degree  Variances  Utilized  in  Computation 
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Figure  5,  Standard  Deviation  of  4v  versus  Revolution 
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Figure  6.  Standard  Deviation  of  A  fa  versus  Revolution 
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Figure  7.  Standard  Deviation  of  Ai  versus  Revolution 
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Figure  9.  Standard  Deviation  At  versus  Revolution 


effects  of  these  terms  on  the  orbit.  Therefore,  the  curves  obtained  by 

summing  from  n  =  7  to  50  correspond  to  the  assumption  of  a  perfectly  known 

gravity  model  through  the  6th  degree  terms  whose  effects  have  been  included 

in  the  orbit  prediction  and  those  obtained  by  summing  from  n  =  13  to  50 

correspond  to  the  assumption  of  a  perfectly  known  gravity  model  through  the 

12th  degree  terms.  Since  the  degree  variances  have  been  calculated  from 

free-air  anomalous  gravity  data,  the  n  =  2  and  n  -  4  degree  variances  do  not 

include  the  preponderant  portion  of  the  and  terms  of  the  terrestrial 

gravity  field.  Thus,  the  n  =  2  to  50  case  corresponds  to  an  orbit  prediction 

using  good  approximations  to  J ^  and  J4  and  an  exact  GM.  The  decision  to  set 

the  upper  limit  of  the  summation  to  50,  while  somewhat  arbitrary,  may  be 

justified  from  examination  of  Table  1.  This  table,  which  was  computed  as  sum - 

2 

ing  a  constant  degree  variance  of  10  mgal  for  all  terms,  presents  the  contri¬ 
bution  to  the  variance  of  the  orbital  elements  of  selected  values  of  n,  normalized 
with  respect  to  the  value  associated  with  n  =  2,  at  the  end  of  16  satellite  revolu¬ 
tions.  The  fall-off  of  each  of  the  values  with  increasing  n  appears  to  justify  the 
truncation  of  the  summation  at  n  =  50.  Since  a  satellite  in  a  circular  orbit  at 
an  altitude  of  100  n  mi  will  complete  approximately  sixteen  revolutions  in  one 
day,  Figs.  4  through  9  yield  the  standard  deviation  in  the  prediction  of  the 
orbital  elements  during  an  interval  of  approximately  one  day.  The  uncertain¬ 
ties  in  the  orbital  elements  can  be  translated  into  radial,  along-track,  and  cross¬ 
track  position  uncertainties  by  noting  that  Ap  is  the  radial  component,  Ai  and 
A 12  can  be  transformed  into  the  cross -track  component  for  polar  orbits,  and 
At  can  be  transformed  into  the  along-track  component.  One-sigma  values  of 
these  components  for  a  one-day  orbital  prediction  resulting  from  the  use  of  an 
incomplete  gravity  model,  as  obtained  from  Figs.  4  through  9,  are  sum¬ 
marized  in  Table  2. 
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APPENDIX  A 


DERIVATION  OF  EQUATIONS  FOR  THE  CORRELATION 
OF  THE  FORCING  FUNCTION 


If  the  ith  component  of  gravity  is  specified  on  the  surface  of  a  sphere 
of  radius  rQ  as  gA  (8,  X)  where  0  and  X  are  the  opherical  coordinates, 
colatitude  and  longitude,  respectively,  and  the  value  of  the  ith  component  of 
gravity  located  at  the  point  0'  ,  X'  with  respect  to  the  point  8  ,  XQ  in  the  0,  X 
coordinate  system  is  gi(8'  ,  X' ),  then  the  covariance  between  values  of  the  ith 
component  separated  by  an  angular  distance  i|>,  C(i|/),  is 


ci<«  =  ir  /  f-s  f%\-  VM*  x',ax'  ■“  W,  <A-'> 

(9.X)  ° 

where  the  gravity  covariance  is  defined  as  the  average  of  all  possible 
products  of  the  ith  gravity  components  separated  by  the  arc  distance  v|i. 

The  radial  component  of  the  covariance  is  most  easily-derived  and  will 
be  considered  first.  If  the  gravipotential  expansion  U  is  considered  to  be 


U(r,  0.  X) 


J  P  (cos  0)  cos  m(X  - 
nm  n  '  '  ' 


(A~2) 


where 

p  =  the  product  of  the  universal  gravitational  constant  and  the  mass  of 
the  earth 

a  -  the  earth's  equatorial  radius 
r  =  the  vehicular  separation  from  the  geocenter 


A-i 


-  - 

*$4 


=  Legendre  polynomials  of  the  first  kind  of  degree  n 


Pm  =  Legendre  associated  functions  of  the  first  kind  of  degree  n 


n 


J  and  Jm  -  coefficients  of  zonal  and  non-Tional  harmonics,  respectively, 
n  n 

then  the  radial  component  gy  of  anomalous  gravity  on  the  .ohere  is 


and  order  m 


(  °°  n 

gr(0. M  =  |—  =  X>  +  D(~)  [w=°8  e) 

t  ' n-2  o 


+  W]j 

m=l  ’ 


(A-3) 


with  a  similar  expansion  for  gy  (6' ,  X' ) 

g  ,e'.x')  =  -*i  X>+  ')(f)°k(eo,xo)P,(c,.  e, 

( n=2  o'  L 


+  y*  J'  (0  ,  X  )  P*"  (cos  e')  cos  m(X'  -  X') 
£mt.  nm  o  o  a  nm' 

m=l 

(A-4) 


where  the  expression  for  J'  (6  ,  X  )  will  be  determined  later.  From  Eq. 
r  n  q  o 

(A-l)  the  averaging  process  over  the  X'  coordinate  is  only  performed  on 
g  (6' ,  X').  Hence,  all  non-zonal  terms  will  average  to  0.  Taking  the  summa¬ 


tion  outside  the  integral  yields  for  (A-i]> 


Cr<4>)  = 


£  t,  P»(«*  «  /”.(2  «r<8o-  V V*  de 


0  n-2 


dX 
o  o 


(A-5) 
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The  relationship  /or  can  be  established  from  the  following  considerations. 
At  every  point  on  the  sphere  gr < 0' ,  X.' )  -  gr(6,  X);  therefore 


P"  fj[,  Pn(c08  0/ }  gr(6'  *  k' >  dS'  =  g ff  Pn(cOB  °' >  «r(6*  V)  dS  {A_6) 


where  S  is  the  spherical  surface  and  dS  is  an  element  of  area.  From  the 
orthogonality  properties  of  ^(cos  6},  the  left  side  of  (A~6)  gives 

- Jn  +  i 


Now 


P„(co.  e'l  =  .Fn(co,  0)  Fn(co.  0O)  +  2  £  in  l  mil 

m=l 

•  P™(co«  6)  P™(cos  0  )  cos[rn(X  -  Xq)]  (A-7) 


If  we  let  J  _  cos  mX  -  A  .  J  _  sin  mX  =  B  ,  then  integration  of 

nm  nm  nm  nm  nm  nm  “ 

the  right  side  of  (A-6)  yields 


Jn<“o'  V  -  W“*  So>  *  .£  t 


mi 


P^icoa  6  ) 
nm  n  o 


cos  mX 

o 


4  B  P  (cos  9  )  sin  mX  ] 
nm  n  '  o'  a 


(A-8) 


Define  the  normalization  factor  for  spherical  harmonics  as 


N2(a’m) 


(A-9) 


A- 3 


where  6„rt  =  ]!.  H  i*  the  Kronecker  delta, 

mo  JO  if  m*0 

Then,  substituting  (A-8)  into  (A-5)  allows  the  integration  to  be 
performed  by  inspection  yielding 


A-4 


Eq,  (A- 11)  may  be  expressed  in  terms  of  Kaula's  degree  variances  by 


cxm 


W  s  MW' 


(T  P  (COS  lb) 
n  n'  ■ 


(A- 14) 


The  covariance  of  the  horizontal  components  of  the  gravity  is  con¬ 
siderably  more  difficult  to  derive  by  the  methods  used  to  obtain  the  covari¬ 
ance  of  the  radial  component.  However,  Kaula  (Ref.  3)  gives  the  covariance 
of  the  deflection  of  the  vertical  in  ter  rut'  of  the  component  p  parallel  to  the 
great  circle  arc  joining  the  wo  point®  and  the  component  t  perpendicular 
to  the  line  joining  the  two  points.  These  components,  in  the  case  of  interest 
here,  represent  the  along-track  and  the  cross -track  component,  respectively. 
Noting  that,  to  a  very  good  approximation,  the  horizontal  component  of  the 
anomalous  gravity  field  is  just  the  total  gravity  force  multiplied  by  the 
defection  angles  (in  radians),  Kaula's  results  can  be  used  to  obtain  for  the 
along-track  covariance  at  the  earth's  surface 


w  ,  '  ,,  •>  f  P?(coa  i|jH 

=  2  [p«(co8  ^ "  "^nTT  j 


(A- 15) 


and  for  the  cross -track  covariance 


cwW0  - 


irp  (n  +  l)n 

tz  2 -  4 


l*rx 


■  Pn-i<cos  +) 

P»-l'co“  **  +  mfiTTT)- 


(A- 16) 


2  2 

where  (coa  iJj)  and  Pn_j(cos  ip)  are  Legendre  associated  functions.  Kaula 
obtained  the  above  results  in  terms  of  cr^  by  differentiating  the  expression 
relating  geoid  height  to  the  geopotential  and  squaring.  On  a  sphere  of  radius 

*Q»  the  expression  for  the  radial  covariance  of  degree  n  is  given  by  the 

(*jl  \2n*i  4  1 

and  the  expression  on  the  surface.  Hence,  to  convert 


Kaula'a  results  to  a  sphere  of  arbitrary  radius  each  terra  of  the  summation 
of  the  right  hand  sides  of  Eqs.  (A-15)  and  (A-16)  must  be  multiplied  by 


However,  it  is  noted  (see  e.  g.  ,  Ref.  5)  that 


Pn(cos  40  =  2Z'Zn(Zk)(2n  :  kk)  C08(n  -  2k)  +  (A- 19) 

k=0  ' 

and  that 


2  3^  r(k+  ~)  r(n  -  k  +  4) 

P»"=°-  «.|£  2'  1 - V 

kTo  r  (f)k!  (n  -  k  -  2)! 

cos(n  -  2k  -  2)  ~  cos(n  -  2k  -  4)»{j  -  j  cos(n  -  2kH‘|  (A-20) 

Thus,  one  sees  that  ail  he  covariance  functions  can  be  expressed  as  linear 
combinations  of  cos  (X.40  terms  where  X  is  an  arbitrary  integer.  In  other 
words,  substitution  of  Eqs.  (A- 19)  and  (A-20)  into  Eqs.  (A-14),  (A-17),  and 
(A  -18)  yielde  covariance  components  of  the  form  given  by  Eq.  (6-f)  in  the 
text. 


A- 6 


APPENDIX  B 


EXPRESSIONS  FOR  ELEMENTS  APPEARING  IN  THE  COVARIANCE 
MATRIX  OF  THE  ORBITAL  ELEMENTS 


Lower  bounds  for  the  covariance  of  the  orbital  elements  are  obtained 
from  direct  integration  of  Eq<  (7)  where  the  elements  of  the  integrand  are 
given  by  Eqs.  (6-d),  (6-e),  and  (6-f).  The  resultant  solution  is  given  by  the 
symmetrical  matrix,  Eq.  (8),  in  the  text. 
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O 
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v  sin  1 
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— TT  'z66 ) 
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and 

Xii  = - jr - 5  I2  I  COB  *<A[(4  -  (L  -  A)  sin  A  -  cos  A  ] 

2(X*  -  i)  L  ° 

+  co»  £  -  40)[cot(4  -  4Q)  -  cos  XA  cos (4  - 
+  X  sin  XAjcos  A[ -2(1,  -  4Q  -  A)  +  sin  2(4  -  4 
-  5  sin  A  -  cos(4  -  4Q  -  2A)»in(4  - 
+  sin(4  -  40)cos(4  -  40  -  A)} 

+  2X2{-(4  -  4q  -  A)  sin  A  cos  XA 

+  sin  (4  -  4Q£sin(4  -  4Q)  -  cos  XA  sin(4  - 


+  X3  sin  XAjcos  A[  2(4  -  4Q  -  A)  -  sin  2(4  -  4q 
-  ai»(&  -  &0)cos(4  -  4q  -  A) 

+  cos(4  -  -  2A)sin(4  -  40  -  A)fJ 


-  A)]+  if 

A)J 

^o  “  A) 


C0  *  A)]  | 

-  A)]  +  sin  A 
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[4  (1  -  COD  XA)  +  4\  sin  XA[  (4  -  4Q  -  A)  -  sin (4  -  4q)] 

+  2X  { cos  XA[  (4  -  4Q  -  A)  sin  A  -  sin(4  -  4Q)sin(4  -  4Q  -  A) 
+  4-2  cos (4  -  4Q  -  A)  -  2  cos  A  ] 

+  sinZ(4  -  ;Q)  -  2  +  2  cos(4  -  4o)| 

+  X3  sin  XA[-  2(4  -  4q  -  A) (4  +  cos  A)  -  7  sin  A 

+  sin(4  -  4q  -  A)cos(4  -  40-  2A) 

-  sin(4  -  4q)cos(4  -  4q  -  A)  +  4  sin(4  -  4Q  -  A) 
+  8  sin(4  -  4q)  -  cos  A  sin  2(4  -  4Q  -  A)] 

4 

+  2X  |  cos  XA[  -  (4  -  ;o  -  A)  sin  A  -  2  +  2  cos(4  -  4q  -  A) 

+  cos  A  -  cos (4  -  40)c°s(;  -  4o  -  A)] 

+  1  +  cosZ(4  -  tQ)  -  2  cos(4  -  40)f 

5  . 
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-  2  sin  (f,  -  iQ  -  A)cos  (4  -  4Q  -  2A) 

+  2  sin(4  -  ycos(t,  -  4Q  -  A) 

+  4  cos (4  -  t,0)sin(£,  -  40  -  A)1 
+  4  sin(4  -  40)  -  4  sin  2(4  -  4Q)} 

+  \b  sin  XAl-  3(4  -  40)2  +  3(4  -  4C)[2  sin(4  -  4Q) 

+  2  Bin  (4  -  4o  -  A)  +  A] 

-  6A  sin(4  -  4q)  +  2  cos(4  -  4Q) 

t  L  cos(4  -  4q  -  A)  -  2  cos  A  -  2 

-  8  sin(4  -  4Q)sin(4  -  40A)}] 
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■* 


\ 


- - ■«  fzl  it  t  -  A)  sin  A  cos  XA  +  ain2(t  -  C  ) 

2(X*  -  if  L  •  -o 

-  cos  XA  sin(4  -  4Q)sin(4  -  4Q  -  A)] 

+  X  sin  XAjcos  A  [  -2(4  -  4Q  -  A)  -  sin  2(4  -  4Q  -  A)] 

-  3  sin  A  +  sin(4  -  4Q  -  A)cos(4  -  4Q  -  2A) 

-  sin  (4  -  40)cos(t  -  i0  -  A)} 

2 

+  2X  [  -  (4  -  4Q  -  A)sin  A  cos  XA  -  cos  A  cos  XA  +  1 

+  cos2(4  -  4Q)  -  cos  XA  coa(4  -  4Q)cos(4  -  4q  -  A)] 


+  X3  sin  XAjcos  A[  2(4  -  4q  -  A)  +  sin  2(C,  -  4  -  A)] 
-  sin  A  -  sin(4  -  4Q  -  A)cos(4  -  40  -  2A) 
+  sin  (t  -  t0)co8(t  -  40  -  A)}] 


Y33=X11 

X-,  - - 5 - 1  I2  sin  XA  sin  &  -  &«,> 

36  2x  (x2  -  ir  L 


+  2X  {  cos  XA[-(4  -  ~  A)  sin  A+  sin(4  -  4Q)sm(t  -  4C  -  A) 

+  cos  (4  -  4  -  A)  +  cos  A] 

-  1  -  cos(4  -  4C)  -  «inZ(4  -  4q)} 

+  X2  sin  XAjcos  A [2(4  -  4Q  -  A)  +  sin  2(4  -  4C  “  A)] 

-  2  sin(4  -  4q  -  A)  -  4  sin  (4  -  4q) 

i 

+  5  sin  A  -  sln(4  -  40  “  A)cos(4  -  4Q  -  2A) 

+  sin(4  -  40)‘0«(;  -  -  A)} 
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+  2X3{co8  XA[  (4  -  4Q  -  A) sin  A  +  cos(4  -  40)(cos(&  -  40  -  A) 

-  cos(4  -  40  -  A)] 

+  coa(4  -  40>[  1  -  cos(4  - 

+  X^  sin  XA  |  -  cob  A[2(4  -  4q  -  A)  +  sin  2(4  -  4C*  A$] 

+  sin(4  -  C0  -  A)eos(4  -  40  "  2 A) 

-  siu( (,  -  40)coa(4  -  4C  -  A)  +  2  sin(4  -  4Q) 

+  2  sin(4  -  4Q  -  A)  -  Bin  A{j 

— T  [-  3(1  -  cos  XA)[  1  +  ccs(4  -  4Q)] 

ir  L 

+  3X  Bin  XA[-  (4  -  &Q  -  A)cos(C  -  4Q)  +  •!»(&  -  &QI  +  A] 

+  X2  |  (?,  »  tQ)[  -  3  sin  (4  -  4Q)  +  3  cob  XA  ain(l  -  40  -  A) 

-  4  sin  A  cos  XA] 

+  cos  XA  [  -  6  +  7  cos  A  +  3  cos(4  -  40  -  A) 

-  6  cos(4  '  4p>  +  4A  8in  A 
+  4  cos(4  -  40)cos(4  ~  *»0  "  A)J 

+  3  cob(4  -  4Q)  -  4  co82ii  -  ;o)  -  !} 

+  X3  Bin  XA  |(4  -  40>U  cosfc  -  4U)  <■  3  cos(4  -  tQ  -  A) 

+  4  cos  A] 
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-  6A  cob(£  -  WQ)  -  °  Bxn^  “  *o' 


.  %  /  a  4  A  4  1  ^  HXI\  A 

-  3  Bin(C,  -  l  -  a;  -  - - 


-  2  sin  2(£  -  &Q  -  A)cos  A 

+  2  cosU  -  £0  -  2A)ain(£,  -  4Q  -  A) 

-  2  sin(^  -  tp)co s((  -  t0  -  A)} 

+  X4  {  a  -  6q)[  3  oin(£  -  l0)  -  3  cos  XA  sin(£  -  tQ  -  A) 

+  4  sin  A  cos  XA] 

+  cos  XA  [3-3  cos  A  -  3  cost’i  *  !>0  ”  A) 

+  3  cob(4  -  40)  -  4A  sin  A 
+  4  sin(t  -  t0)uin(4  -  tQ  -  A)] 

-  4  sin2(^  -  C0)[ 

+  X5  sin  XA  J  -  -  £q)[3  cos(C,  -  CQ)  +  3  cos(-  £  -t0-  A) 

+  4  cos  A] 

+  3 A  cos(;  -  C,o)  +  3  sin  {£  -  &0)  -  5  sin  A 
+  3  sin(t  -  -  A)  +  3A  +  4A  cos  A 

+  2  sin  2(t  -  tc  -  A)cos  A 
+  2  sin(r,  -  L  )coo«t  -  4  -  A) 


2  cos  (4  -  ;Q  -  2A)sin(4  -  iQ  -  A)j] 


v  -  »"  * 
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fciit 


i 

I 


+  cos  coi  4  ■  cos  XA  cos(4  -  A)]} 

+  sin  \A  j  cos  A  [  -  2(4  -  4Q  "  A)  -  sin  2(4Q  +  A) 

+  sin  2(4  -  A)]  -  4  sin  A 

-  sin  40  cos(4Q  +  A) 

+sin(4c  +  A)cos(40  +  2A)  +  sin  4  cos  (4  -  A) 

-  ain(4  -  A)cos(4  -  2A)J 


+  2X2|-  (4  -  CQ  -  A)  sin  A  cos  XA 

+  sin  4J  sin  4Q  -  cos  XA  sin(4Q  +  A)] 
+  sin  4(sin(4  -  cos  XA  sin(4  -  A)]} 


+  X3  sin  XA  {  cos  A[2(4  -  4Q  -  A)  +  sin  2(4Q  +  A) 

-  sin  2(4  -  A)]  +  sin  4Q  cos(4Q  +  A) 

-  sin(40  +  A)cos(40  +  2A)  -  sin  4  cos(4  -  A) 
+  sin(4  -  A)cos(4  -  2A)|j 
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[  -2  sin  cos  -  2  sin  t  cos  C 

y  w  u 

+  cos  XA  {  3  sin  5q  coa(£o  +  A)  +  3  sin  £  cos(£  -  A) 

-  *in(£0  +  A)coa(£o  4  2A)  -  sin(5  -  A)cos(5  -  2A) 
+  2  sin  A[  sin2(Z,  -  A)  -  sin2(r,0  +  A)]} 

4  X  ain  XA  {  -  sin  A[  sin  2(t,Q  +  A)  +  sin  2(5,  -  A)] 

4  2  COB  A[sin2(C0  +  A)  -  sin2(C,  -  A)]} 

4  X2|^2  sin5,a  cos  £0  +  2  sin  £  cos  C 

4  cos  XA  {  -  sin  cos(5q  +  A)  ••  sin  5  cos(£  -  A) 

4  ain(5Q  4  A)[cos(C0  +  2A)  -  2  cos  CQ] 

4  sin(C  -  A)[cos(5  -  2A)  -  2  cos  £] 

4  2  sin  A[  sin2(5.Q  +  A)  -  sin2(C  -  All] 

4  X3  sin  XAjsin  Afsin  2(£q  4  A)  4  sin  2(5  -  A)] 

4  2  cos  A[  sin2(C  -  A)  -  sin2(;Q  4  A)]f) 


j^2(4  -  4Q  -  A)sin  A  cos  XA 

+  2  siii  4Q[  gin  4Q  -  cos  XA  sin(4Q  +  A)] 

+  2  sin  4[  sin  4  ”  cos  XA  sin  (4  -  A)] 

+  X.  sin  XAj  cos  A[  -  2(4  -  4Q  -  A)  -  sin  2(4  -  A) 

+  sin  2(4q  +  A)]  -  4  sin  A 

+  sin40  C03(40  +  A)  -  Bin  4  cos(4  -  A) 

-  sin  (4q  +  A)cos(40  +  2A)  +  sin(4  -  A)coa(4  -  2A)| 

+  X^  {  -  2(4  -  4Q  -  A)  sin  A  cos  XA 

+  2  cos  4  fcos  4^  -  cos  XA  cos  ( l ,  +  A)] 
o  o  y 

+  2  coo  4[cos  4  “  cos  XA  cos(4  -  A)]( 

+  X3  sin  XAjcos  A[2(4  -  V>0  “  A)  +  sin  2(4  -  A)  -  sin  2(40  +  A)] 

-  sin  4q  cos(4q  +  A)  +  sin  4  cos(4  -  A) 

+  sin(4Q  +  A)cos(40  +  2A)  -  sin(4  -  A)cos(4  -  2A)[| 


R-ifJ 


\*  -  VUB  IV4%/  -  AM  A41  »AU 


4  3  X2{3(4  -  4Q)2(i  -  coa  XA) 

4  cos  \A[3A2  4  4-8  coa(4  -  40)J  -  4 
+  8  a-i(&  -  40)f 

4  3X3  sin  XA[  2(4  -  4Q)3  -  3A(4  -  4Q)2 

4  8(4  -  CD  -  A)co$(4  -  r,o)  4  A3  4  4A 

-  8  sin(4  -  C0)l 

4  2X4{-  9(4  -  i0)2(1  "  cos  XA> 

4  (4  -  4o)[  12  sin(4  -  4U>  -  12  coa  XA  sin(4  -  4( 
4  8  sin  A  cos  XA] 

4  coa  XA[-  9A2  4  i5  -  20  coa  A  -  12  cos(4  -  4( 
4  24  coa(4  -  40)  -  8A  3i«  A 
-  8  coa(4  -  40)c°*U  -  40  -  A>1  +  5 

*y 

-  12  cos (4  -  t0)  +  8  cos  (4  -  40)| 

4  2X5  sin  XA!-  6(  4  -  40)3  +  9A(4  -  4Q>2 

-  4(4  -  40)[  6  cos(4  -  4q)  4  3  cos(4  -  4, 

4  2  coa  A] 


.  3A3  +  15A  +  24A  cos(4  -  40>  +  24  sin^  “  to) 
.32  sin  A  +  12  ain<4  -  40  -  A)  +  8A  cos  A 
+  4  sin  2(5,  -  ;o  -  A)cos  A 

-  4  cos  (C  -  4C  -  2A)sin(4  -  40  -  A) 

+  4  sin(4  -  4q)cob(4  -  4Q  -  A)| 

+  V6  {  9(4  ■  4q)2(  1  -  cos  XA) 

+  (4  -  40)[  -  24  sin{  4  -  4C)  +  24  cos  XA  sin(  4  -  t,Q-  A) 
-  16  sin  A  cos  XA] 

+  cos  XApA2  -  24  +  24  cos  A  +  24  cos(4  -  4Q  “  A) 

-  24  cos(4  -  4q)  +  i6A  sin  A 

-  16  sin(4  -  40)*in(4  -  tQ  -  A)] 

+  16  sinZ(4  -  40){ 

+  X7  sin  XA{6(4  -  4a)3  -  9A(4  ••  4e)2 

+  8(4  -  40)[J  cos (4  -  40)+  3  co»(4  -  4q  -  A) 

+  2  cos  A  ] 

+  3A2  -  24A  cos (4  -  40)  -  24  sin(4  -  4Q) 

+  32  sin  A  -  24  sin(4  -  4Q  -  A)  -  24A 

-  16A  cos  A  -  8  sin  2(4  -  4Q  -  A)cos  A 
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-  8  #in(!i  -  £c)cos(t  -  ^  ~  A) 

+  8  cos (4  -  -  2A)sin(&  -  t,Q  -  A)jj 

Z66  =  Z55 

The  indeterminacies  occuring  in  the  above  equations,  when  X.  -  0  and  when 
X  =  1,  are  resolved  by  direct  application  of  L'Hospital's  rule. 
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